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ABSTRACT 

We propose a statistical method to search for fluctuations in the temperature of the intergalactic 
medium (IGM) using the Lyman a forest. The power on small scales (~ 25 km/s) is used as a ther- 
mometer and fluctuations of this power are constrained. The method is illustrated using Q1422+231. 
We see no evidence of temperature fluctuations. We show that in a model with two temperatures that 
occupy comparable fractions of the spectra, the ratio of small scale powers is constrained to be smaller 
than 3.5 (corresponding to a factor of 2.5 in temperature). We show that approximately ten quasars are 
needed constrain factors of two fluctuations in small scale power power. 

Subject headings: cosmic microwave background — methods: data analysis 
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1. INTRODUCTION 

The Lyman a forest has become one of the major tools 
of cosmology. It has been used to set constraints on many 
of the parameters of the cosmological model, such as the 
amplitude and slope of the spectrum of primor dial density 
fluctuations and the nat u re of the dark matter ( Croft et al] 
19981 ICroft et al. 1999|, |Whitc fc Croft 2000|, Narayanan! 



et al. 20001 , |Zaldarriaga, Hui fc Tegmark 2000| , |Ci : oft eg 
2000 ^ In the popular cosmologies, the forest a rises 



ratlicir^ iaturally and both cosmolog i cal uimulatioii 



et al. 



1994 [Hcrnquist et al, 1995|, [Zhang et al, 1-9*05 



trate it using the spectrum of one quasar and investigate 
what type of constraints on the temperature fluctuations 
could be obtained with more data. 

A possible explanation for the high temperature of the 
IGM around redshift z ~ 3 is heat deposited in the gas 
by the double ionization of He. Depending upon the na- 
ture of the source responsible for the reionization of the 
universe and its dumpiness, Hell reionization can be sub- 
stantially delayed re lative to H reionization and could hap- 
pen around z ~ 3 ( Miralda-Escude fc Rees 1994 , Madau 
& Meiksin 1994 Miralda-Escude et al. 20*00|}T 



Miralda-Escude et al. 1996|, Muecket et al. 1996] , |Wadsley 
fc Bond 10961, f heuns et al. 1*198) as well as analytical 



models have been used to understand its properti es 
et al. 1992|, jRciscncggcr fc Miralda-Escude 1991) 



Th a atatiatica 



T 

by th e physical properties of the intergalactic medium 



Davidson 19971, IGncidn fc Hui 1996| , pott ct al. 1997| , Hui 
fc Gncdin 19974 [Hui et al. 1997b| ) 



A natural consequence of this scenario are fluctuations 
in the temperature of the IGM around z ~ 3. In fact 
studies of the Hell Lyman a forest show that the mean ab- 
sorption is increasing rapidly with redshift for redshifts be- 
tween 2 and 3. Moreover several gaps of size ~ lOOOkm/s 
where the transmitted flux is high have been observed 



propcrtica of the forest flux arc affected 



( Reimers et al. 1997, [Anderson et al. 
Bmcttc ct al. 



2000 



2000)7 



1998, Heap et al 



(IGM) such as its temperature or the slope of its equation 
of state affect. In turn these physical properties of the IGM 
are sensitive to the ionization history of the universe and 
the cosmological parameters. An intensive effort has been 
devoted to understanding what different scenarios predict 
for the forest and in const raining those scenarios with the 
available observations (eg. Hui fc Gncdin 1997a[ Gncdin fc 
Hui 1998|, [Bchaye et al. 1999| , [McDonald ct al. 200q , [Croft 



et al 



Auniform IGM temperature has been assumed in most 
of the theoretical and observational Lyman a forest work 
so far. We live in an inhomogeneous universe, inhomoge- 
neous heating of the IGM due to a late Hell reionization 
is only one of the possible sources of temperature fluctu- 
ations. It is timely to investigate possible ways of detect- 
ing temperature fluctuations. In this paper we use the 
small scale powe r spectrum of the for e st flux as a measure 
of temperature (Thcuns ct al. 2000, Zaldarriaga, Hui & 



Th d temperature of the IGM has been determined in 



Tegmark 2000) an introduce the tools necessary to search 



two different ways, by measuring the wid ths of the lines 
dSchaye et al 199S| , pVIcDonald et al. 2000| ) and b y study- 



ing the small scale power sp ectrum of the flux (Zaldar- 
riaga, Hui & Tegmark 200C). Both methods find val- 
ues of the temperature around T ~ 20, 00 K which are 



higher than theoretical expectations (ie. Hui fc Giiedin| 
19974 Hui et al. 1997b[ ). Both these type of stud les as 
well as most theoretical models assume that the tempera- 
ture of the IGM is constant in space. The aim of this paper 
is to formulate a method for testing this assumption, illus- 



for fluctuations in the small scale power. 

The paper is organized as follows, in section §|2j we sum- 
marize the properties of the small scale power spectrum 
and its dependence on the IGM temperature, in §|3| we 
introduce the tools we use in the search for temperature 
fluctuations, in jfl we apply our technique to Q1422+231 
|Kim ct al. 19971), a quasar at redshift z = 3.62 and in 



we estimate how much data is needed to substantially 
improve our constraints. We conclude in §H 
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TrangmigQmn power spectrum from Zaldarriaga 
Tegmark 2000. Models with have varying l'o [Tfj = 
650)(km/s) 2 increases from top to bottom, 7 was set 



(150, 250, 

to 0.2 and the mean transmission was set to 0.7). The data points 
are from McDonald et al. (1999). 



2. THE POWER SPECTRUM AS A MEASURE OF THE IGM 
TEMPERATURE 

The power spectrum of the flux in the Lyman a for- 
est can be used as an indication of the IGM temperatur e 
flThcuns ct al. 2000| , [Zaldarriaga, Hui fc Tegmark 20001 ). 
Figure [l] shows the predictions of a set of models with 
different temperatures together with the measurements of 
McDonald et al. 19991 The simulations used in this p aper 



are described in Zaldarriaga, Hui fc Tegmark 20001 We 



use PM simulations to solve for the dark matter and then 
use analytic scaling relations to construct mock flux spec- 
tra. Our simulations have 128 3 particles in a box that at 
z — 3 is 3200km/s in size. 

Given that the small scale power is sensitive to the tem- 
perature of the IGM it is reasonable to try to detect fluc- 
tuations in the IGM temperature by searching for fluc- 
tuations in the small scale power. To do so we have to 
measure the power in different positions along the spec- 
trum and/or on lines of sight towards different quasars. 
We could use Fourier analysis on pieces of the spectra but 
we prefer to use wavelets, a natural tool to obtain mea- 
sures of power that are simultaneously localized in real 
and frequency space. In the next section we review the 
properties of the wavelet expansion that are relevant for 
our work. 

3. WAVELETS 

Wavelets have been used in a number of statistical stud- 



ies of the forest. In Mciksin 200C wavelets were introduced 
for data compression and several statistical properties of 
the wavelet coefficients were presented. 



For our purposes th e more relevant work is that of The- 
mis & Zaroubi 2000. Just as we will do in this paper 



the authors used wavelet coefficients to measure the small 
scale power in the spectrum. They presented the cumula- 
tive probability distribution of the squares of the wavelet 
coefficients and focus most of their effort on comparing 
models with different temperatures and equations of state 
and thus different small scale power. They conclude that 
the wavelet coefficients in these different models have dif- 
ferent enough distributions to be able to tell the mod- 
els ap art observationally. In Zaldarriaga, Hui & Tegmark 



2000 we used Fourier analysis and measurements of the 



small scale power in the literatu re to constrain the equa- 
tion of state of the IGM. Both iThcuns & Zaroubi 2000 



and Zaldarriaga, Hui fc Tegmark 200q are very similar in 
this respect, they propose the same "thermometer" for the 
IGM, the small scale power in the forest. 

Furthermore, Thcuns & Zaroubi 200C argue that be- 
cause wavelet coefficients are localized in real space they 
can be used to measure fluctuations in the temperature. 
They showed results for the cumulative distribution of a 
model with two temperatures and compare them to models 
with one single temperature. Unfortunately it is unclear 
how much of the difference is due to the models having a 
different mean temperatures and how much is due to the 
fluctuations themselves, or how much data is needed to get 
tight constraints on fluctuations. Our work in this paper 
buil ds on the ideas presented in Thcuns fc Zaroubi 2000 



and Zaldarriaga, Hui & Tegmark 2000 with special focus 
on statistical tests of temperature fluctuations. 

In this section we introduce the concepts necessary to 
understand the use of wavelets in the context of this paper. 

3.1. Definitions 

The best way to introduce the necessary concepts is by 
looking at the wavelet decomposition of a piece of spec- 
tra from Q1422+231 flKim et al. 1997[ ). In figure | we 
show 3200 km/s of the spectra together with wavelet coef- 
ficients. This is approximately a 10th of the Lyman a part 
of the spectra. We will use all the spectra later, here we 
chose to show a small portion mainly for clarity but also 
because this is the length of the simulations we will have 
to compare with. The wavelet transform of the spectra 
is obtained by convolving the spectra with a set of test 
functions which we will call ip n (x — Xi), 



a n {xi) = / dxF(x)ip n (x - Xi), 



(1) 



x denotes position along the spectra and F(x) is the flux. 
The test functions are characterized by two numbers, the 
index n and Xi. Increasing values of n correspond to 
higher spatial frequencies, Xi determines the position of 
the wavelet along the spectrum. The values of de- 
termined by the wavelet scheme. For the wavelets we are 
using there are 2™ values of Xi in each order n. Figure || 
shows a few examples of these wavelets to illustrate how as 
n increases, the wavelet corresponds to higher spatial fre- 
quencies. In this work we use the Daubenchie s 4 wavelets 



imple mented in a Numerical Recipes routine (Press ct al 
1992D. 



For our purposes the wavelets are just a set of filters 
that we apply to the data. They have the advantage of 
being localized in both Fourier and spatial domain, of be- 
ing an orthogonal and complete set, and perhaps most 
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Fig. 2. — Wavelet transform of 3200km/s of the spectrum of 
Q1422+231. On the top we show the actual spectrum. In the 
middle panels we plot the amplitude of the wavelet coefficients along 
the spectra for four different orders. In the bottom panel we show 
examples of the wavelet functions if) for the orders whose coefficients 
are shown. 
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Fig. 3. — In the top panel we show the power spectrum of the 
forest in a cosmological simulation measured both using Fourier 
and wavelet analysis. The line corresponds to Fourier analysis, the 
points are the equivalent measures using wavelets. On the bottom 
we show the window functions for the n = 4, 5, 6, 7 wavelets. The 
power measured using wavelets is a weighted average of the Fourier 
power spectrum with the windows as weights. 



importantly that packages that perform the transform are 
readily available. 

To further de-mystify the wavelets we can look at the 
power spectrum of the forest from a wavelet perspective. 
Rather than computing squares of Fourier coefficient we 
can obtain a measure of power by averaging the squares of 
the wavelets coefficients of a fixed n but different positions 
along the spectrum. We have 



1 N 



(2) 



the constant of proportionality depends on the normaliza- 
tion convention of ip n . The wavelets are well localized in 
real space so they must have a broad response in Fourier 
space. These response or window functions for a few orders 

are shown in figure H Ry construction window functions of 



succe ssive n differ only by a factor of 2 scaling in the x axis. 
In figure [| we also show the power spectrum measured in 
the conventional Fourier way or using wavelts. The x co- 
ordinate of each wavelet measurement is set to the peak 
of the window function, but the windows are clearly very 
broad. 

We should point out that the n label is arbitrary in the 
sense that the label for a wavelet that probes a particular 
range of spatial scales depends on the length of the spec- 
trum being analyzed. For example, if we start with a piece 
of spectrum that is twice as long (6400 km/s) all the labels 
would be shifted by one. That is, the new n = 7 wavelet 
probes the same scales as the old n = 6 did. As long as 
one is consistent with the normalization convention the 
new wavelet coefficients of order n + 1 will be equal to the 
old coefficients of order n at the same position along the 



spectrum. To avoid confusion we will use the labels cor- 
responding to 3200km/s in all our discussions even when 
we analyze the full length of Q1422+231 spectrum. Figure 
|H indicates exactly what scales we probe when we discuss 
each order n. 

It is useful to discuss how the power in the wavelet 
coefficients depends on the temperature. In figure || we 
illustrate the rate of change of power with temperature 
obtained in our simulations. Clearly the smallest spatial 
scales are most sensitive. The figure also shows that for 
n = 7 a factor of 2 change in the temperature, which 
could be expected in models were He has been doubly 
ionized in patches of the spectrum, could lead to a fac- 
tor of 2.5 change in the observed power. For the figure 
me assumed an isothermal model, the temperature of the 
gas is independ ent of its density. In Zaldarriaga, Hui & 



Tegmark 2000 we showed that the small scale power is 



sensitive to the temperature over a range of overdensities 
0.9 < p/p < 1.8. The ratio of temperatures in figure || 
should be interpreted as the ratio of temperatures on that 
range of overdensities. This fact can be particularly im- 
portant in the case of energy injection by He reionization 
as the process will not only change the temperature of the 
gas but also its equation of state, so the relevant param- 
eter is by how much the temperature is increased in the 
range of overdensities that are probed by the small scale 
power. 

In summary wavelets of different scales n measure power 
on different spatial frequencies. Their frequency response 
is broad, the trade-off for the wavelets being localized in 
real space. The wavelets probing the small scales are very 
sensitive to the temperature of the IGM and we will use 
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Fig. 4. — Ratio of power (average of the square of the wavelet 
coefficients) for models with different temperatures normalized to a 
model with T = 300(fcr?i/s) 2 = 1.8 X 10 4 K in our simulations. We 
show the results for n = 6 and n = 7 wavelets. 



these coefficients as probes of the temperature. We will 
take advantage of the localization in real space of the 
wavelets to search for temperature fluctuations. 

3.2. Statistical properties of the wavelet coefficients 

Now that we have defined the wavelet coefficients we 
will turn to some analysis of their statistical properties in 
the case of the Lyman a forest flux. 

In the rest of the paper we will be interested in the small 
scalepower, that we will be using as a thermometer. In fig- 
ure H we show a histogram of the amplitude of the n = 7 
coefficient obtained in our simulations. The variance of 
this distribution proportional to the flux power spectrum 
(equation ||). In the the top panel we show the histograms 
for two different temperatures, T = 400 (800) (km/s) 2 = 
2.4 (4.8) x 10 4 K. Clearly the model with the largest tem- 
perature has the smaller small scale power, the variance of 
the distribution is smaller. 

In the bottom panel of figure || we have rescaled the x 
axis by dividing it by the variance. The two histograms 
now fall on top of each other, although there are some 
small differences. For comparison we also show the his- 
tograms of the coefficients of Q1422+231 which have a 
distribution which is remarkably similar to the coefficients 
in the simulation. 

In figure we show a Gaussian together with the 
histograms of the coefficients in the simulations and in 
Q1422+231. There are clear differences between the Gaus- 
sian and the rest. The amplitudes of the wavelets of the 
flux have a probability of being near zero which is much 
larger than one would expect in a Gaussian of the same 
variance. This is clearly noticeable in the example of fig- 
ure ||, several stretches of the spectrum can be identified 
where the amplitude of the wavelet coefficients are very 
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Fig. 5. — In the top panel we show the histograms of the ampli- 
tude of n = 7 wavelet coefficients in our 3200km/s simulations for 
two different temperature, T = 400(km/s) 2 (larger variance) and 
T = 800(km/s) 2 (smaller variance). In the bottom panel we have 
rescaled the axis using the variance. We also show the distribu- 
tion for the coefficients of Q1422 + 231 (dotted line) and a simple 
Gaussian (solid line). 



small. The distribution of the coefficient on the simulation 
also has tails, it is much probable to get large coefficients 
than in a Gaussian distribution. 

To lay the foundation for our method to probe for tem- 
perature fluctuations we also need to study the spatial cor- 
relations of the wavelet coefficients. In figure ^ we show the 
correlation coefficient of the wavelet amplitudes (CC) as 
a function of separation along the spectra (Sv) for several 
scale indices n, 



CC(Sv) = 



< a n (x)a n (x + Sv) > 



<ai> 



(3) 



The correlation coefficient decays fast with separation and 
then levels off at a value below 1%. In our future anal- 
ysis we will be considering only the smaller scales, corre- 
sponding to n — 6 or 7. The correlations in these cases 
become negligible (below 1%) for separations larger than 
- lOOkm/s. 

3.3. Analytic tools 

In order to detect temperature fluctuations we want 
to look for variations in the statistical properties of the 
wavelet coefficients as we look at different parts of the 
spectra or as we compare different lines of sight. It is clear 
that any statements about fluctuations will be statistical 
in nature, the variance in the distribution of the wavelet 
coefficients is our thermometer, not the coefficients them- 
selves. Thus we cannot measure the temperature on a 
pixel by pixel basis, we need to average over many pixels. 

Our measure of power is the square of the wavelet co- 
efficient, thus we will introduce averages of the square of 
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Fig. 6. — Correlation coefficient for wavelets of different indices 
n as a function of separation Sv. For the smallest scales (n < 6 the 
correlations become negligible (< 1%) foe separations larger that 
~ lOOkm/s. 

these coefficients which we will call w, 

P i 

The sum is done over contiguous portions of the spectra 
of length L. Different portions do not overlap, and there 
are N p coefficient in each chuck of size L. N p depends on 
the order n one is considering, N p — 2 n (L/3200km/s). By 
definition the mean of w is the measure of power and thus 
temperature. The mean only depends on the index n and 
not on L. In figure |we show the histogram or probability 
distribution (PDF) of w for several choices of L in our 
simulation. We show the histograms for two temperatures. 
It is clear from the figure that the model with the larger 
temperature has a smaller mean w. 

In figure || we have rescaled the x axis in figure [?] using 
the mean of w (id) for each temperature. The histograms 
for the two temperatures are quite similar, we are thus en- 
couraged to find some simple analytic model. Our statistic 
w is an average of squares of many variables which as fig- 
ure ^| shows quickly become uncorrelated. We could expect 
that a x 2 distribution would be a good model for the ob- 
served distributions. The effective number of degrees of 
freedom v of the x 2 distributions need to be determined. 
Furthermore because the different panels correspond to 
values of L that increase by factors of 2, the effective num- 
ber of degrees of freedom should increase by factors of 2. 
In figure we plot this simple analytic approximation, 

P(w\w) = - V e- 
I (v/2)2w 

u = wv/2w, (5) 

where for the different panels we choose v — 1,2,4,8. 
Changing slightly the values of v improves or worsens the 



Fig. 2_ — Probability distribution function of w (defined in 
equation H) for four different choices of L and n = 7. Each 
panel shows two temperatures, T = 400(km/s) 2 (larger mean) and 
T = 800(km/s) 2 (smaller mean). The mean of w is independent of 
L, it depends on the temperature and n. 

agreements depending on the panel. We choose these val- 
ues doing a simple x by eye of the above plots. Our con- 
straints on the temperature fluctuations are insensitive to 
the details of this choice, different values of v give the same 
answer, as long as the analytic formulas fit the panels rea- 
sonably well. Note also that the distribution for the data 
of Q1422+231 is also quite similar, there are some differ- 
ences for the largest smoothing window but in that case 
there are only a few samples Q1422+231, so the discrep- 
ancy is not significant. 

We will use this analytic form for the pdf of w through- 
out the rest of our work. Although the approximation is 
not perfect our analytic formulas give a good description 
what is observed in the simulations. We are going to be 
searching for variations in this distribution as we move 
along the spectra, so it is very useful that we have an 
analytic form and furthermore that this analytic form de- 
pends on only one parameter (w). When formulated this 
way, the question of fluctuations in the temperature can 
be very simply stated. Does a model with a single value of 
w fit the data better than a model were w varies in some 
particular way across the spectrum. 

We could take a different approaches to look for fluctu- 
ations in the statistical properties of the coefficient. We 
could directly compare the observed distributions with 
what is obtained in simulations where a uniform tempera- 
ture was assumed and search for differences. The problem 
with this approach is that one must establish that this dif- 
ferences are not a consequence of systematic effects intro- 
duced in the simulation process by some missing physics 
or some numerical artifact. We believe that our method is 
more robust as it does not use simulations directly, only 
to motivate our choice of distribution. We are essentially 
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Fig. 8. — Rescaled probability distribution function of w for four 
different choices of L and n = 7. Each panel shows two temper- 
atures, T = 400(km/s) 2 (solid) and T = 800(km/s) 2 (dashed) to- 
gether with our analytic model. In addition in the L = 400km/s 
panel (the smoothing scale we use in our statistical study) we show 
the pdf for Q1422+231. 



looking for changes in the parameters of the distribution 
we use to model the data as we go along the spectrum. 
Any change in this parameters is a signal of fluctuations 
even if our distribution is not a perfect model of the data. 

4. CONSTRAINTS ON TEMPERATURE FLUCTUATIONS 

Our objective is to use the small scale power as a ther- 
mometer, and this thermometer to search for fluctuations. 
In the previous sections we have established a model for 
the probability distribution of the square of the wavelets 
coefficients. This model is dependent on one parameter, 
which we chose to be the mean power on that scale. The 
objective of this section is to describe the formalism nec- 
essary to constrain spatial fluctuations in the mean power, 
such as those that would be produced by fluctuation in the 
temperature. 

Let as call w the vector with the square of the ampli- 
tude of the wavelet coefficient as a function of position 
along the spectra. Let us denote N the dimension of that 
vector, which is given by the ratio of the length of the spec- 
trum to length of the averaging filter in equation For 
simplicity we will call each of the entries of these vectors 
pixels, although they are a combination of several pixels 
in the initial spectrum. 

In the previous sections we found a model for P(wi\wi), 
the probability of measuring Wi given that mean Wi of w 
at position i along the spectrum, (equation]^). To describe 
our model for the fluctuation we need to introduce another 
vector (w) that specifies the value of id along the spectrum. 
The dimension of w is also N. The probability of measuring 



P(w) = / dwP(w)P(w|w) 



JV 



dwi ■ ■ ■ dw n P(w) JJP(wi|t0i). 



(6) 



i=l 



Different models for the temperature distribution differ in 
what P(w) is. For example for a uniform temperature 
model where for every pixel the mean power is w we have, 

P(w) = 5 D ( Wl -w)---6 D (w N -w) (7) 

By combining equations (jjj), (||) and (^) we obtain the 
probability of measuring w under the assumption that the 
temperature is uniform (P"), 

lnP" = ^2hxP(wi\w) 



i 

N\n[2wT{v/2)/v], 



(8) 



where Ui — Wiv/2w. 

We want to constrain fluctuations in the temperature 
by comparing the probability of having measured a vector 
w in a particular model for the temperature fluctuations 
(P f ) with the same probability in a uniform temperature 
model (P™). The ratio of and P" is the likelihood ra- 
tio for the two models. In order to compute P^ we need 
to specify the value of P(w), our model for the tempera- 
ture fluctuations. Clearly there are many different choices 
of P(w) which describe different physical scenarios and it 
would be impossible for us to try every possibility. Instead 
in this section we have chosen some particular examples 
and found constraints in those scenarios. 

In all our models we will assume that there are only two 
different temperature values possible and that the mean 
power for each of these two temperatures wa and wb, 
differ by a factor j3 = wb /wa- We will also define / as the 
fraction of the spectrum with mean wa ■ The mean power 
averaged over all pixels becomes, 



w = fw A + (1 - f)w B - 



(9) 



When comparing different models we will keep this mean, 
w fixed to the observed value in the spectrum. We are not 
interested in comparing models with different mean power 
but models with or without spatial temperature fluctua- 
tions in the power. 

There are two probability ratios that define the prob- 
lem, 



In 



P(Wj\w A ) 

P(wi\w) 



-«<(l-/)(l-/3) 



In' 



\w B ) 



P(wi\w) 



-Uif 



/?(!-/)) 
/ + /?(! 



/3 



(10) 



where Ui = WiV /2w. The expressions for P-^/P" in the dif- 
ferent models for the fluctuations will always be expressed 
in terms of these two quantities. 



4.1. Independent pixels 

We first consider a model in which every pixel is inde- 
pendent of the rest and can be in one of two temperatures, 
corresponding to mean powers wa and wb- The prior on 
the small scale power becomes, 

P(w) - \\[f5 D { Wl - + - f) 5 °^i - (11) 
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The likelihood ratio p//P u is, 



ln[PVP"]=5Jln[/ 

i 

+ (i-/) 



P(Wj\w A ) 
P(Wi\w) 

P{Wj\w B ) -, 
P(Wi\w) 1 



This expression remains constant if we transform (3 
and / - (1 - /). 



(12) 



1/(3 



We use equation (12) and the full spectra of Q1422+231 
to put simultaneous constraints on / and (3 which are 
shown in figure [| We do not need to consider values of 
(3 <1 because they can be obtained from (3 > 1 using the 
symmetry property of the likelihood ratio. We considered 
models with a likelihood ratio smaller than 0.05 as ruled 
out. As can be expected the tightest constraints on (3 are 
obtained when / ~ 0.5. With the data we have we can 
say that (3 < 3.5 (/ - 0.5). 

As could have been guessed the tightest constraints are 
obtained when both temperatures occupy a significant 
fraction of the spectrum, / ~ 0.5. However it is interesting 
to note.that the contour is not symmetric around / = 0.5. 
Figure g shows that it is harder to put constraints when / 
is near zero than when it is near one. The reason for this is 
that when / ~ and (3 > 1, most of the spectrum has the 
larger small scale power wb (low temperature) and only a 
small fraction has the lower power wa (high temperature). 
Even when w is high, figure || shows that it is more likely 
that in any given pixel the value of w is low, thus it is 
hard to detect the small excess of pixels with low ws that 
is produced by the small fraction of the spectrum that has 
w A - 

4.2. Sharp transition 

We now consider a case in which there is a sharp tran- 
sition between wa and wb which occurs at pixel number 
itr- The prior becomes, 



p(w) = n * 



W A 



)LK( 



Wi—WB) 



(13) 



%>% tr 



The likelihood ratio P^/P" for this case becomes, 



ln[PVP"]= £ ln[ P{Wi ^ A) 

i<i tT 

+ E M 



l>l tT 



P{w t \w) 

P(Wj\w B ) , 

P{w t \w) 1 



(14) 



For this model we define / = itr/N . We compute /P" 
for all possible values of it r and plot the constraints on (3 
in figure M. We plot the results in terms of / rather than 




Fig. 9. — Two dimensional constraints on the fraction / of spec- 
trum with mean wa and ratio /3 = uig/wA- The dark regions have 
a small likelihood ratio, P-^/P" < 0.05. The panels correspond to 
different models for the fluctuation, independent pixels (top left), 
sharp transition (top right), fixed size bubbles of size 1200km/s 
(lower left) and of size 2400km/s (lower right). 



itr- For this model we plot constraints from both (3 bigger 
and smaller than one, as both cases are not equivalent. 
One corresponds to a jump from low to high temperature 
and the other to the opposite case. 

These model is the one where the tightest constraints 
can be obtained because it is the less random. Given i tr 
the temperature of all pixels is specified. We can can con- 
strain 0.5 < (3 < 2 for almost all values of /. 

4.3. Fixed size bubbles 

Finally we will consider a case in which there are bub- 
bles of a fixed size with power wb and the rest of the 
spectrum has wa- To compute the P(w) under this sce- 
nario we reinterpret equation (^|) as an average over P(w) 
of P(w|w). We can use a Montecarlo technique to perform 
this average. We create realizations simply by laying down 
Nb bubbles of a fixed size in the spectrum. We choose ran- 
domly the centers of the bubbles. For each realization we 



compute the likelihood using equation (10) and the aver- 
age that over all realizations to compute the P(w). The 
model has two free parameters, the size of each individual 
bubble and the number of bubbles (or equivalently /, the 
fraction of the spectrum having wa and wb)- 

In figure || we show the constraints on / and (3 for two 
particular bubble sizes, 1200 km/s and 2400 km/s. In this 
model (3s larger and smaller than one are not equivalent. 
Our results for independent pixels can be interpreted as a 
constraint when the typical size of the bubble is 400 km/s. 
As can be seen in figure ^ our constraints end up being very 
similar for all bubble sizes. For / ~ 0.5 the ratio of powers 
have to be in the range 0.3 < (3 < 3.5. 

5. FUTURE 

In the previous sections we have developed a method for 
constraining temperature fluctuations. We illustrated our 
method using data from Q1422+231. The constraints on 
the ratio of the powers obtained are somewhat larger than 
what one might expect to be present in the universe. In 
this section we compute the number of spectra needed to 
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substantially improve the constraints. 

The statistical question we will ask is what constraints 
on (3 can be obtained from a certain amount of data in the 
independent pixel model under the assumption that the 
universe has uniform temperature. Equation (12) implies 
that the ratio of the probabilities of a uniform model to 
the fluctuating temperature model (P-^/P 11 ) is, 



ln(PVP u )=E ln [/ 

i 

+ (i - /) 



P{w l \w A ) 
P{wi\w) 

P(Wi\w B )-, 



In [/ 



P(Wi\w) 

N J dwP(w\w) 
P(w\wa) 



P(w\w) 



(1-/) 



P(w\w B ) 
P(w\w) 



(15) 



In the last line we approximated the sum by a mean over 
the distribution from which the Wi are assumed to be 
drawn, P(w\w), and we used their independence. The last 
line of equation (15) is a measure of distance between two 
probabili ty distribution functions called Kullback-Leibler 
entropy ( |Kullback 1959| ). 

We can set the ratio in equation ( |l5| ) to 0.05 and solve 
for N as a function of (3 and /. Weshow the results in 
figure 10. We plot only points for j3 > 1, other values 
of (3 can be obtained using the transformation (3 — > 1/(3, 
/ — > (1 — /). We have changed variables form TV to N q , 
the number of quasars with the same number of pixels as 
Q1422+231. The figure shows that ten quasars are needed 
to obtain constraints on (3 ~ 2 for / ~ 0.5. The figure 
also illustrates the fact that it is very difficult to obtain 
constraints that are much tighter than a factor 2 in (3. 

It is interesting to understand why it gets so much 
harder to obtain constraints when (3 approches one. Fig- 
ure [l^ shows that the number of quasars scales as N q oc 
l/(/3 — l) 4 . The key point is that we are comparing mod- 
els that by construction have the same mean w (the same 
mean temperature). Thus what we are actually compar- 
ing is the variance of the w distribution in both models. 
We can compute how many samples are needed to distin- 
guish two models that have different variance. We define 
the variance as v —< w 2 > —w 2 , its value depends on the 
model (fluctuating or uniform temperature), but w is the 
same across models. We are assuming all the pixels are 
independant so we can estimate the variance as, 



-E 



W: 



w 



(16) 



The question becomes how many samples (N p ) are 
needed to distinguish a model with uniform temperature 
with one characterized by (3 and /. We define the signal 
to noise ratio as, 



(S/N) 2 = [v f -v u ] 2 /vai(v), 



(17) 



the ratio of the difference in variance between the model 
with fluctuations and the one without to the variance of 
the variance (which we calculate in the uniform model). 
In our example the signal to noise ratio becomes, 
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Fig. 10. — Number of quasars similar to Q1422+231 needed to 
obtain a constraint for several values of / as a function of /3 = 
wb/wa- 



{S/Nf 



cN, 



(l-/3) 4 / 2 (l-/) 2 
[/ + 13(1 - /)] 4 



(18) 



where the constant c depends on the shape of the distri- 
bution of w and is defined as c =< w 2 > 2 /(< w 4 > — < 
w 2 > 2 ). The scaling of the signal to noise with the number 
of pixels is (1 — /3) 4 as figure || shows. It is the fact that 
we are looking at the difference in the variance of w for 
models constructed to have the same w, what makes dis- 
tinguishing the models so difficult as (3 — > 1. This fact is 
what forces N q to scale as 1/(13— l) 4 instead of l/(/3 — l) 2 . 

Finally with more data we could not only get better con- 
straints on model parameters such as (3 but also constrain 
more ambitious models. For example, rather than have 
fixed size bubbles we could draw the bubble sizes from 
some distribution or allow for correlations between their 
positions. Moreover with more data one has to study in 
detail certain sources of systematic error in our method 
like contamination by metal lines or the effect of evolution 
along the spectra. 

6. CONCLUSIONS 

The uniformity of the temperature in the IGM is an 
untested assumption in most theoretical modelling of the 
Lyman a forest. We have presented a statistical method 
to test this assumption. Our method uses the small scale 
power of the forest flux as a thermometer and searches 
for fluctuations in this power. As a result the method 
is statistical in nature. An important advantage of the 
method we suggest is that it does not require a direct 
comparison with numerical simulations. It searches for 
changes in the statistical properties in different parts of 
the spectra. 

We illustrated our method with Q1422+231. We as- 
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sumed that in different regions of the spectrum the IGM 
could be in one of two temperatures. We found that the 
ratio of small scale power in these two type regions could 
not differ by more than a factor of 3.5 if the two type of re- 
gions occupy a comparable fraction of the spectrum. This 
ratio of power translates in a factor of approximately 2.5 
in temperature. If the transition between the two regimes 
is sharp, the constraints are more stringent, the ratio of 
powers cannot differ by more than a factor of 2. 

We have calculated how much data is needed to improve 
this constraints substantially. We showed that ten quasars 
are needed to constrain fluctuations of factors of two in 
the power, or equivalently factors of 1.7 in temperature. 
Furthermore with more data, more ambitious models for 
the fluctuations can be considered. 

Acknowledgements The author acknowledges very 
useful comments and suggestions from Joop Schaye and 
John Bahcall. Support for this work is provided by the 
Hubble Fellowship HF-01 116-01-98 A from STScI, oper- 
ated by AURA, Inc. under NASA contract NAS5-26555. 



Wadsley J. W., Bond J. R. 1996, Proceeding of the 12th Kings 



astro.ph/9^12148 



Conference, eds. Clarke D., West M., RASP, 
White M., Croft R. A. C. 2000, preprint ^stro-ph/UUU124Y| 
Zhang Y., Anninos P., Norman M. L. .19957 
Zaldarriaga M., Hui L., Tegmark M., istro-ph/0011559 



r5T 



REFERENCES 

Anderson et al. 1998, preprint a,stro-ph/9808105| 
Bi H. C, Boerner C, Chu Y.,LW2, A&A, 2bb, 1 
Bi H. C, Davidsen A. F. 1997, ApJ, 479, 523 
Cen R., Miralda-Escude J., Ostriker J. P., Rauch M 
L9 

Croft R. A. C, Weinberg D. H., 
488, 532 

Croft R. A. C, Weinberg D. H., 
495, 44 

Croft R. A. O, Weinberg D. H., 

1999, ApJ, 520, 1 
Croft R. A. C, Weinberg D. H., 

Katz N., Kirman D., Tytler D 
Gncdin N. Y., Hui L., 1996 
Gnedin N. Y., Hui L., 1998, MNRAS, 296, 44 
Heap et al. 2000, ApJ, 534, 69 

Hcrnquist L., Katz N., Weinberg D. H., Miralda-Escude J. 1995 

ApJ, 457, L5 
Hui L., Gnedin N. Y. 1997a, MNRAS, 292, 27 
Hui L., Gnedin N. Y., Zhang Y. 1997b, ApJ, 486, 599 
Kim T. S., Hu E. M., Cowie L. L., Songaila, A. 1997, AJ, 114, 1 
Kullback S., Information Theory and Statistics, Wiley, New York 

1959 

Madau P., Meiksin A., 1994 ApJ433 L53 
McDonald P., Miralda-Escude J., Rauch M 
Barlow T, ,a. , Cen R., Ostriker J. P., 



1994, ApJ, 437, 
Hernquist L., Katz N. 1997, ApJ, 
Katz N., Hernquist L. 1998, ApJ, 
Pettini M., Hernquist L., Katz N. 
R^lto M , Rurlsc S, i Hernquist L., 



astro-ph/0012324 
ApJ, 472, /3 



ph/c 911196 
McDonald FT 
Barlow T 



., Sargent W 
1999 



preprint |astro-| 



Miralda-Escude J., Rauch M. 
i. , Cen R., Ostriker J. P., 



W. 



ph/( 005553 . 

Meiksin A., 2000 preprint astro-ph/0002148 
Miralda-Escude J., Rees M. iyy4, MiNRAS 



442. 



, Sargent 
2000 preprint astro 



266, 343 



Miralda-Escude J 
582 

Miralda-Escude, J., Haehnelt, M. & Rees 
Muecket J. P., Petitjean P., Kates R. E., 
308, 17 

Narayanan V K , Spergel D. N., Dave R., 

astro-ph/0001247 
Press W. H., leukolsky S. A., Vetterling 



Cen R., Ostriker J. P., Rauch M. 1996, ApJ, 471, 



M. J. 2000, ApJ, 530, 1 
Ricdiger R. 1996, A&A, 

Ma C. P. 2000, preprint 



_ W. T., Flannery B. P., 
1992, Numerical Recipes in Fortran (2ed; New York: Cambridge 



Univ. Press) 
Reimers et al. 1997, AA, 327, 890 

Reisenegger A., Miralda-Escude J. 1995, ApJ, 449, 476 
Schaye J., Theuns T., Leonard A., Efstathiou G. 1999, MNRAS, 
57 



310, 



ill;, 



■ a m 



fripp T. M., Jenkins E. B.. 



Smette A., Heap S. R., 

Songaila A., preprint |astro-ph/0012193| 
Theuns T., Leonard A., Eistathiou G., Pearce F. R., Thomas P. A. 

1998, MNRAS, 301, 478 
Theuns T., Schaye J., Haehnelt M. 2000, MNRAS 
Theuns T., Zaroubi S., preprint ^stro-ph/0002172 



315, 600 



